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We develop one numerical method to compute black and gray solitons with Kerr-type nonlo- 
cal nonlinearity. As two examples of nonlocal cases, the gray soliton with exponentially decaying 
nonlocal response or with Gaussian nonlocal response are discussed. For such two nonlocal cases, 
the analytical form of the tails of nonlocal gray soliton is presented and the analytical relationship 
for the maximal transverse velocity of nonlocal gray soliton to the characteristic nonlocal length is 
obtained. 
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I. INTRODUCTION 

In present years spatial solitons with Kerr-type non- 
local nonlinearity have attracted a great amount of 
studies. It is indicated the nonlocality of the nonlin- 
earity ensures the existence of stable multidimensional 
solitons [U, 0, S B]- The nonlocal bright spatial soli- 
tons have been experimentally observed in nematic liquid 
crystaljl, 0, Q and in the lead glassjl]. The propaga- 
tion and interaction properties of nonlocal bright/black 
solitons are greatly different from that of local soli- 
tons. For example, the dependent functions of the beam 
power and phase constant on the beam width for non- 
local bright solitons are very different from those for lo- 
cal bright solitons d, [13, El; There exists higher order 
Hermite- Gaussian like (1+1) dimensional bright nonlo- 
cal solitonsjl, [13, [l2| and Laguerre-Gaussian like (1+2) 
dimensional bright nonlocal solitons [ol. [ill. [T^ ITsl. IT4. ITU : 
Bright nonlocal solitons with tt phase difference attract 
rather repel each other [l^, [H, [l3i 13; Nonlocal black 
solitons can attract each other and two black solitons 
can form a bound state [H, H^, [2l[. The nonlocal gray 
solitons also can form a bound state [2^. 

In this paper we present a numerical method to com- 
pute black and gray soliton solutions with Kerr-type 
nonlocal nonlinearity. This numerical method is a gen- 
eralized version of the so called spectral renormaliza- 
tion method presented by M. J. Ablowitz and Z. H. 
Musslimani[22| with which they compute the bright soli- 
ton solutions. Our numerical method can be easily gen- 
eralized to compute a large family of gray soliton so- 
lutions with other type of nonlinearity. It is worth to 
note that Yaroslav V. Kartashov and Lluis Torner[2^ 
has rather successfully studied the nonlocal gray soliton 
with an exponential decaying nonlocal response in a dif- 
ferent manner. They has revealed that the gray soliton 
velocity depends on the nonlocality degree and pointed 
out the maximal velocity of gray soliton monotonically 
decreases with the degree of nonlocality [24]. In our pa- 
per, besides of the nonlocal case of an exponential de- 



caying nonlocal response, the nonlocal case of Gaussian 
nonlocal response is also discussed. The analytical form 
of the tails of nonlocal gray soliton is presented. It is 
indicated that nonlocal gray solitons can have exponen- 
tially decaying tails or have exponentially decaying oscil- 
latory tails. The analytical relationship for the maximal 
transverse velocity of nonlocal gray soliton to the charac- 
teristic nonlocal length is presented. It is indicated when 
the characteristic nonlocal length is less than some criti- 
cal value, such maximal transverse velocity is a constant 
and equal to that of local gray soliton, otherwise such 
maximal velocity will decrease with the increasing of the 
characteristic nonlocal length. That mostly agrees with 
the results of reference [24] with slight but not trivial ex- 
ception. It was said in reference [24] that such maximal 
velocity monotonically decreases with the characteristic 
nonlocal length. But as will be shown in our paper the 
maximal velocity does not vary with the characteristic 
nonlocal length when the characteristic nonlocal length 
is less than some critical value. Our paper poses a prob- 
lem for such inconsistence. 



II. NUMERICAL METHOD TO COMPUTE 
NONLOCAL GRAY SOLITONS 



The propagation of an (1+1) dimensional optical beam 
in Kerr-type nonlocal self-defocusing media can be de- 
scribed by this following (1+1) dimensional nonlocal non- 



linear Schodinge r eq uati on( NNLSE) [1| 
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where u{x, z) is the complex amplitude envelop of the 
light beam, \u{x,z)\^ is the light intensity, x and z are 
transverse and longitude coordinates respectively, R{x), 
{J R{x)dx — 1) is the real symmetric nonlocal response 
function, and n{x,z) — — J R{x — C)l'"(^i ^)I^'^C is the 
light-induced perturbed refractive index. Note that not 
stated otherwise all integrals in this paper will extend 
over the whole x axis. When R{x) — 6{x), equa- 
tion ([T]) will reduce to the local nonlinear Schodinger 



2 



equation(NLSE) 



and 



.du 1 d'^u , ,9 

which has black and gray sohton solutions [25 

u{x, z) = V(x)e*['^^+'^(^)1 



11/2 



(f){x) ~ r\\/ 1 — B'^x + arctan 



B tanh(77Sa;) 



(2) 



(3a) 
(3b) 

(3c) 
(3d) 



It is easy to prove that i\j{x) ^^+°°> ^ 4)' [x) ^~^+°°> 

ryVl — -B^ = So we have /3 = — (^ +^7^) which, as 
will be shown, also applies to nonlocal gray soliton. 

In this paper we numerically compute the black and 
gray soliton solutions of NNLSE ([T]) that take these fol- 
lowing form 

m(x, z) = V(a;)e*''^+''^("') (4a) 
?/)*(a;) = ■(/'(a;) lim i\)(x) = rj > (4b) 

[3* = (3 (/)*(x)=0(x) (4c) 

Substitution of Eq. ([4a|) into ([T]), we have an equation for 
■0(2;) and 0(a;). Since both the real and imaginary part 
of such resulting equation must vanish, we have 

- /3V' + \r - ^(V-OV - V- / R{x - Oi^'im = (5) 



2(/)V + V-'V = 



(6) 



From Eq. ([51), we have -^{tlj^cj)') = 0, which results in 
ip^(t)' — const. We set this constant as iirj"^, where ji* = ^ 
is another constant. So we have 



(7) 



Since il^{x) °°> 77 as assumed in Eq. (j4bp . from 

Eq. (O, we have (f)' {x) ^ — which in turn results 

in 4){x) c + fix, where c is a constant. The phase 

jump through the soliton is defined as 2c. Substituting 
Eq. dll) into (O, we get 

-(3tp + 1^" - ^ - V; 1 i?(x - c)^'(e)de = 0. (8) 



Since "0(2^) 



X — ^+00 



rj and J R{x)dx = 1, from Eq. 



when X +00, we have —(3r] — — 77^ = that leads 



to 



(9) 
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/ i?(a;-O0'(eK = O. 



(10) 



So the gray soliton takes a form of u{x, z) = 

0(x) exp[— 7(-^ + rf')z + i(j){x)], where 0(x), 4>{x) satisfy 
Eqs. (fTO|) and It is easy to prove that the GaHlean 
transformed solution 

ui{x,z) = M(a; + ^z, z)e*^^^^^^^^^^ 

= 0(x + ^2)e'[-'''^+"^(^+''^)"^("'+^^)l 

= 0l(a; + Ai^)e-*'''^ (11) 

where 0i(a;+/xz) = ^(a;+^2)e»[</'(^+/'2)-/'(^+M^)]^ also sat- 
isfies the NNLSE HI). Since |Mi(a;,z)p = |0i(x-fAi^)|2 = 
02(x-|-/iz), the gray soliton moves in velocity — /i with re- 
spect to the coordinate system. On the other hand, since 
0(x) ^^+°°> and 0(x) ^^^°°> c + /ix, from Eq. (|11|). 

we have ui(a;,z) '^"*~^°°) jye'^'^"'' ^•'j which implies the 
background of gray soliton is at rest with respect to the 
coordinate system. So the soliton ui{x, z) moves in trans- 
verse velocity —fj, with respect to the background inten- 
sity of gray soliton. In this paper we simply refer —^l as 
the transverse velocity of the soliton. 

Assuming 0(x) °°) 77', from Eq. ()10|) . we have rj' = 

2 4 

±77. When ji ^ 0,the term requires ■0(a;) ^ 0, which 
in turn requires 77' = 77. So the case n ^ corresponds to 
the gray solitons and 0(— a;) = 0(a;), where we have set 
a; = be the symmetric center. 

We study the case 7^ in this section and leave the 
case /i = to be discussed in the next section. Let 



0(a;) = r]-x{x), 



(12) 



where xi^x) = x{x) and x{x) ^^^°°> Q. Then Eq. ([TO] 
turns into 



2 ^ 2^ ^ 2 



iv - x)^ 

-iv-x)J R{x ~ Oix'iO - ^rjxiOn = 0.(13) 

We discrete the function x(a;) in Xj = x(^^+(i^l)'^2;)i 
where —h < x < h is the sample window, Ax is the 
sample step and 1 < j < n = + 1 . Define the discrete 
Fourier transform (DFT) ^ by 



X, =^[xL-=EFjkXk, (14a) 

n 

X. ='^"Mxb- = EFlkXk, (14b) 
fe=i 



3 



where Fjk = ^exp[i^(j — l)(k — 1)]. Performing the 
DFT on Eq. we have 



+ -fx, + 



(?7 - 



where 17,- 



we have 



iv-x) Rix-OixHO-^r^xiOn 



Let 



0, 

(15) 



|^ 2sm[|0-l)] y 



2 2 ^ 2 



= A0(a;), (16) 



(?? - Ae')3 



= 0. 
(17) 

Projecting Eq. pT| onto 6*, we obtain an equation for A 
^e(A) 



.7 = 1 



3A?76'2 - 37j^e - A' 



2/J31 



(?7 - A6i)3 



= 0. 



iv^X9) / i?(a;-e)[A0'(O-2770(e)]rfe 



(18) 



On the other hand, from Eq. ^T7} . we get 



2 



^ [(77 - xe) J R{x - Oi^o^iO - ^voion] . 



2 



(19) 



where r is a positive constant. 

We use Eqs. (fT8|) and (fTg]) to iteratively compute gray 
soHton solutions. For an initial function 0\(x)., e.g. a 
Gaussian function, from Eq. (fT5|) . we get Ai which sat- 
isfies the equation ^^^(Ai) = 0. Then from Eq. we 
get function 6*2 = For 7ti > 1, wc get the iter- 

ation scheme s:^Q^^(Xrn) — 0,6'm+i — !^\^[9m]- Perform 
the iteration until the convergence is achieved. Then we 



jet tp{x) = T] — X9{x) and (j)[x) = 



As an exarnple, we consider this following nonlocal 

case [H, i, i i, i, m, S M [O, M M, M, M 



n — w 



, drn 



(20) 
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FIG. 1: Amplitudes of gray soliton solutions with exponen- 
tial nonlocal response (a) or with Gaussian nonlocal response 
(b). The parameters used are r; = = 1/3 and black line 
corresponds to local case w = red line ui = 1, green line 
It) = 3 and blue line u; = 5 both for (a) and (b). 




FIG. 2: Numerical simulation of the gray soliton with Gaus- 
sian nonlocal response and with parameters 7;; = 5,77 = 1,^ = 
1/3. 
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FIG. 3: Phase jump through the nonlocal gray soliton with 
an exponential nonlocal response. Black line is the local case 
TO = 0; Red line w = 1\ Green line 10 = 3; Blue line w = 5. 
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FIG. 4: (a) Black line is the amplitude of gray soliton so- 
lution with exponential nonlocal response and with parame- 
ters ri = 1, fj. = 0.43, w — 5, and the red line is the fitting 
function 1 - 0.137 exp(-0.0753|a;|) cos(0.596|a;| + 0.055); (b) 
Gray soliton solution with Gaussian nonlocal response with 
parameters 77 = 1,^ = 0.4, ui = 5 and the fitting function 
1 - 0.238 exp(-0.0733la;|) cos(0.713|2;l + 0.107). 



which results in n{x) — — J ^ exp(- 
where R{x) = ^ exp(- 



is the exponential decaying 
nonlocal response function. In another example, we con- 



sider the Gaussian nonlocal response [1, H, H, [13, [O, [3, 
m, [13, Ei i?(a;) = ^ exp • For such two non- 



local cases, w is referred as the characteristic nonlocal 
length. In Fig. ([1]) some nonlocal gray soliton solutions 
are shown. The numerical simulation result in Fig. ([2]) 
indicates the numerical nonlocal gray soliton solution 
obtained can describe the soliton state very well. The 
phase jump through the nonlocal gray soliton is shown in 
Fig. ([3]). As shown by Fig. ([4]) the nonlocal gray solitons 
can have exponentially decaying oscillatory tails for both 
the exponential nonlocal response case and the Gaussian 
nonlocal response case. 

Now we investigate the form of the decaying tails of 
gray solitons. In case |x(a;)| <C 77, to the leading order 
Eq. can be linearized to 

- - Ix" + 2ry' / - OxiOd^ = 0. (21) 



Since Eq. pT|) is linear for x(x), the superposition theo- 
rem applies. 

For the exponential nonlocal response R{x) — 
exp(— — ), Eq. (PTjl reduces to two coupled equations 



1 

2w 



-2m'x - Ix" + 2v^f{x) = 0, (22a) 
fix) - w^f'ix) = x{^)- (22b) 



The solution of (1221) '^^'^ be assumed as 

X(a;) = aexp(Ax). (23) 
Substituting ^ into (QH), we have 



X =4 



1 - A2w2 



Ai' I X 



(24) 



The eigenvalue problem of the above equation provides 
an equation for A 



A^ = 4 



1 - A2u;2 



(25) 



which results in 



2«;2 



(26) 



When w ^ 0, from Eq. ((26l) to the leading order we have 

A^ = 4(77^ - fj,^){l + irj^w^) 
The roots of Eq. ^ are 

A = ±(Ai +^2), 



(27) 
(28) 



where 



Ai 



A2 



4^2 



(29) 
(30) 



For w = 5, r; = = 0.43, we have Ai = 0.0753, A2 = 
0.596 which are used by the fitting function in Fig. ([4]). 

To obtain exponentially decaying tails, A must be a 
real number and we have 



1 - 4w^^^ -f 4^^772 - ^2 > 0, 



-1 + 47i;^^^ + 4^1/772 - /i2 < 0, 
which in turn result into 



< < 77^ for w < — 

477 

4u'77 — 1 , , ^ 1 1 

On the other hand when 



1 - Aw'^i? + Awy/rf - 7^2 > 0, 



-l + Aw^fi"^ +4w\/Tf^^ > 0, 



(31a) 
(31b) 

(32a) 
(32b) 



(33a) 
(33b) 



A is a complex number and exponentially decaying os- 
cillatory tails are obtained. From inequalities l[33|) . we 
have 



n / 2 / 4u;77 -1,2 



for 



w>i-. (34) 



5 




FIG. 5: The area under the red hne is the parameter space 
for gray sohton with exponential decaying nonlocal response 
when 77 = 1. The left area of the blue line corresponds to 
gray solitons with exponentially decaying tails; the right area 
corresponds to gray solitons with exponentially decaying os- 
cillatory tails. 



FIG. 6: The area under the red line is the parameter space 
for gray soliton with Gaussian nonlocal response when 77 = 1. 
The left area of the blue line corresponds to gray solitons 
with exponentially decaying tails; the right area corresponds 
to gray solitons with exponentially decaying oscillatory tails. 



where 



In conclusion, when w < 77- the maximal transverse 
velocity does not vary with the characteristic nonlocal 
length w and is equal to 77 the maximal velocity of lo- 
cal gray soliton; When w > the maximal transverse 

velocity is \J^^^^- The maximal transverse velocity 
will decrease when w > ^ with the increasing of the 
nonlocal length w. Such results are shown in Fig. ([S]). 
It is worth to note that there is slight but not trivial 
inconsistence between our results and the results of ref- 
erence[24]. As shown in figure 4(a) of reference[24] and 
pointed out in reference [24], the maximal velocity mono- 
tonically decreases with the nonlocality degree (the char- 
acteristic nonlocal length in this paper) . But as has been 
pointed out by our paper and shown in figure ([5]) in this 
paper, the maximal velocity does not vary with the char- 
acteristic nonlocal length w when w < 

On the form of the tail of the nonlocal gray soliton 
we arrive at, for w < the gray solitons always have 
exponentially decaying tails; for u; > the gray solitons 
always have exponentially decaying oscillatory tails; for 
^ < w < the gray solitons can have exponentially 

decaying oscillatory tails when < /i^ < '^^^i^ or have 



exponentially decaying tails when 

Now we consider the nonlocal case with a Gaussian 
nonlocal response R{x) — exp(— ^). Here we only 
discuss the exponentially decaying oscillatory tails. How- 
ever the exponentially decaying tails can be discussed in 
a similar way. For the exponentially decaying oscillatory 
tails 



x{x) — exp(— Ax) cos(Ka;), 



(35) 



A = 



B 



cos ■ 



— sm 

K 2 



1 ^^{x^-^"^) w^Xk 

—e 4 sm . 

K 2 



(37) 
(38) 



So Eq. (mi) turns into 

- ix" + ^V^Bx' + {2ifA - 2t?)x = 0. (39) 

The eigenvalue problem of the above equation provides 
two coupled equations for A and k 

A = -2tj^B, (40a) 
K ^ V-2(2772A-2/i2)- 477452 (40b) 

For example, when w = 5, 77 = 1, /i = 0.4, from (|40p we 
get A ~ 0.0733, k — 0.713 which are used by the fitting 
function in Fig. ([4]). 

In a similar way of deducing the parameter space in 
the exponentially decaying nonlocal case, we also get the 
parameter space in the Gaussian nonlocal case. When 
w < 1/77, the maximal transverse velocity is a constant 
equal to rj; when w > l/?7, the maximal velocity is equal 

to ^ i+M'?""'" !. On the other hand when w > I/77, the 
gray soliton always has exponentially decaying oscilla- 
tory tail; When w < l/y/erj, the gray soliton always has 
exponentially decaying tail; when l/y^rj < w < I/77, 
the gray soliton can have exponentially decaying tail for 

77 > /i > a / ^+^"('^ ™ ) or exponentially decaying oscilla- 



we have 



tory tail for ^ < \/ _ These results are shown 

in Fig. dSl). 



III. NUMERICAL METHOD TO COMPUTE 
NONLOCAL BLACK SOLITON SOLUTIONS 



r°° 1 

f{x) = / — ^cxp 
J-00 w^n 

= Ax{x)+Bx'{x), 



(36) 



In the case — Eq. (jTU]) reduces to 

V/' 



^2^+ / i?(a:-Orm = 0. (41) 
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On the other hand, from Eq. p7)) we have 



r9, + r^^e, - [e J R{x - o^pmdi] . 



FIG. 7: Amplitudes of nonlocal black soliton solutions with 
Gaussian nonlocal response. The parameters used are rj = 1 
and black line w = 0, red line w — 1, green line w — 3 and 
blue line w = 5. 



We consider the nonlocal black sohton solutions with 
properties 

ipi-x) ^ -V(a;) ^{x) ±77 ?/>(0) = (42) 



(49) 



Let 



X(^) = ^ i.e. V'(x)= rx(t)dt (43) 
ax Jo 



So xi~x) = Xix) and xi^) 
Eq. (Hll), we obtain 

1 



0. Acting on 



-2i; [ R{x - OmxiOd^ = (44) 

(45) 

-2AVi J R{x - OMmm - 0, (46) 

Where ipi{x) = d{t)dt. Performing DFT on Eq. (jie]) . 
we have 



Let 

x{x) = \e{x), 

and substitute it into Eq. (|44|) . we have 



R{x-0i^U0d^ 



2^-1 / R{x-OMmOd^ 



0. (47) 



Projecting Eq. (|T7)) onto 0, we obtain an equation for A 
i/e(A) 



-^6*, - A^^ 
2 ^ 



2^-1 / i?(x-0^i(0^(Orf^ 



= 



(48) 



The iteration scheme is ^^^(Am) — 0, (?„i+i = 
^x,n[Qm]- Some black nonlocal soliton solutions with 
Gaussian nonlocal response are shown in Fig. ([7|). 



IV. CONCLUSION 



One numerical method to compute black and gray soli- 
ton solutions with Kerr-type nonlocal nonlinearity is de- 
veloped. As two examples of nonlocal cases, the gray 
soliton with exponentially decaying nonlocal response or 
with Gaussian nonlocal response are discussed. The an- 
alytical form of the tails of nonlocal gray soliton is pre- 
sented. It is indicated that nonlocal gray solitons can 
have exponentially decaying tails or have exponentially 
decaying oscillatory tails. The analytical relationship for 
the maximal transverse velocity of nonlocal gray soliton 
to the characteristic nonlocal length is presented. It is 
indicated when the characteristic nonlocal length is less 
than some critical value, such maximal transverse veloc- 
ity is a constant and equal to that of local gray soliton, 
otherwise such maximal velocity will decrease with the 
increasing of the characteristic nonlocal length. 
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